UNCLASSIFIED 

AD 


DEFENSE  DOCUMENTATION  CENTER 

FOR 

SCIENTIFIC  AND  TECHNICAL  INFORMATION 

CAMERON  STATION.  ALEXANDRIA.  VIRGINIA 


UNCLASSIFIED 


NOTICE:  When  government  or  other  drawings,  speci¬ 
fications  or  other  data  are  used  for  any  purpose 
other  than  in  connection  with  a  definitely  related 
government  procurement  operation,  the  U.  S. 
Government  thereby  incurs  no  responsibility,  nor  any 
obligation  whatsoever;  and  the  fact  that  the  Govern¬ 
ment  may  have  formulated,  furnished,  or  in  any  way 
supplied  the  said  drawings,  specifications,  or  other 
data  is  not  to  be  regarded  by  implication  or  other¬ 
wise  as  in  any  manner  licensing  the  holder  or  any 
other  person  or  corporation,  or  conveying  any  rightB 
or  permission  to  manufacture,  use  or  sell  any 
patented  invention  that  may  in  any  way  be  related 
thereto. 


00 

co 

TjH 

*  CO 


ESD-TDR-  63-477 


TM-3641 


A  BRIEF  SURVEY  OF  SATELLITE  ORBIT  COMPUTATION 
METHODS  AND  MAJOR  PERTURBATIVE  EFFECTS 


CD 

CD 


TECHNICAL  DOCUMENTARY  REPORT  NO.  ESD-TDR- 63-477 


t-LJ 


L«J 

CrJ 

o 

_ ! 

TSt 


oo 

-aC 


MARCH  1964 
R.  H.  Greene 


Prepared  for 


496 L  SYSTEM  PROGRAM  OFFICE 
ELECTRONIC  SYSTEMS  DIVISION 
AIR  FORCE  SYSTEMS  COMMAND 
UNITED  STATES  AIR  FORCE 
L.  G.  Hanscom  Field,  Bedford,  Massachusetts 


Project  496 L 
Prepared  by 

THE  MITRE  CORPORATION 
Bedford,  Massachusetts 
Contract  AF19(628)-2390 


Copies  available  at  Office  of  Technical  Services, 
Department  of  Commerce. 


Qualified  requesters  may  obtain  copies  from  DDC. 
Orders  will  be  expedited  If  placedthrough  the  librarian 
or  other  person  designated  to  request  documents 
from  DDC. 


When  US  Government  drawings,  specifications,  or 
other  data  are  used  for  any  purpose  other  than  a 
definitely  related  government  procurement  operation, 
the  government  thereby  Incurs  no  responsibility 
nor  any  obligation  whatsoever;  and  the  fact  that  the 
government  may  have  formulated,  furnished,  or  In 
any  way  supplied  the  said  drawings,  specifications, 
or  other  data  Is  not  to  be  regarded  by  Implication 
or  otherwise,  as  In  any  manner  licensing  the  holder 
or  any  other  person  or  corporation,  or  conveying 
any  rights  or  permission  to  manufacture,  use,  or  sell 
any  patented  Invention  that  may  In  any  way  be  related 
thereto. 


Do  not  return  this  copy.  Retain  or  destroy. 


ESD-TDR-  63-477 


TM-3641 


A  BRIEF  SURVEY  OF  SATELLITE  ORBIT  COMPUTATION 
METHODS  AND  MAJOR  PERTURBATIVE  EFFECTS 


TECHNICAL  DOCUMENTARY  REPORT  NO.  ESD-TDR- 63-477 

MARCH  1964 
R.  H.  Greene 

Prepared  for 

496 L  SYSTEM  PROGRAM  OFFICE 
ELECTRONIC  SYSTEMS  DIVISION 
AIR  FORCE  SYSTEMS  COMMAND 
UNITED  STATES  AIR  FORCE 
L.  G.  Hanscom  Field,  Bedford,  Massachusetts 


Project  496L 
Prepared  by 


THE  MITRE  CORPORATION 
Bedford,  Massachusetts 
Contract  AF19(628)-2390 


TM-3G41 


ABSTRACT 


A  brief  descriptive  survey  of  the  various  satellite  orbit 
computation  methods  is  presented  here  in  order  to  give  the  systems 
engineer  an  introduction  to  the  various  methods  and  problems  which 
currently  exist. 
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A  BRIEF  SURVEY  OF  SATELLITE  ORBIT  COMPUTATION 
METHODS  AND  MAJOR  PERTURBATIVE  EFFECTS 


INTRODUCTION 


SECTION  I 


Many  systems  engineers  who  have  very  limited  knowledge  of  celestial 
mechanics,  but  who  are,  nevertheless,  occasionally  concerned  with  designing 
space  systems,  may  find  a  brief  introduction  to  the  practical  problems  of  com¬ 
puting  satellite  orbits  useful  for  their  work.  The  following  exposition  makes  no 
attempt  at  completeness,  but  merely  discusses  some  of  the  more  common  and 
useful  computation  methods.  The  mathematics  involved  is  kept  to  a  bare  mini¬ 
mum,  since  it  would  be  almost  impossible  to  give  even  a  fairly  complete  listing 
of  the  various  formulae  used  in  these  orbit  computation  methods  without  filling 
several  books.  However,  references  and  bibliography  compiled  at  the  end  of 
this  paper  will  serve  to  guide  those  more  interested  in  details  to  some  of  the 
most  important  sources. 
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SECTION  II 

REFERENCE  SYSTEMS  z 


The  most  commonly  used  reference  system  for  earth  satellites  is  the 
geocentric  equatorial  system  which  is  fixed  in  inertial  space  (see  Fig.  1).  The 
x-axis  points  toward  the  vernal  equinox  y  which  is  that  point  in  inertial  space 
where  the  apparent  orbit  of  the  sun  crosses  the  equator  from  south  to  north, 
the  z-axis  points  toward  the  North  Pole,  and  the  y-axis  forms  a  right-hand 
system.  The  x-y  plane,  which  is  the  equatorial  plane,  is  the  fundamental  plane 
of  the  system. 

The  orbit  of  a  satellite  is  an  ellipse  with  one  focus  at  the  center  of  the 
earth  (see  Fig.  2),  and  is  completely  described  by  six  independent  parameters. 
Three  of  these  parameters  f2  ,  w  ,  and  i  fix  the  orientation  of  the  orbit  in 
space,  and  the  other  three  parameters  a,  e,  and  T  fix  the  size  and  shape 
of  the  ellipse  and  the  time  of  the  point  of  closest  approach  of  the  satellite. 

Using  Figs.  1  and  2  for  illustration,  we  have: 
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a  =  the  semi-major  axis  of  the  satellite  orbit, 
e  =  the  eccentricity  of  the  ellipse, 
i  =  the  inclination  of  the  orbit  plane. 

=  the  longitude  of  the  node. 
uj  =  the  argument  of  perigee. 

T  =  the  time  of  perigee  passage. 


Fig.  2  Reference  System 

o  =  center  of  earth  B  =  apogee  height 

A  =  perigee  height  q  =  perigee  distance 

Normally,  the  element  M,  the  mean  anomaly,  is  used  in  place  of  the  time  of 
perigee  passage. 

M  =  k  a'3/2(t  -T) 
e 

where  ke  is  the  Gaussian  constant. 
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SECTION  III 

THE  DETERMINATION  OF  SATELLITE  ORBITS 

The  determination  of  a  satellite  orbit  may  be  outlined  as  follows: 

(a)  The  first  determination  of  "intermediate  elements"  E  from 

1 

observations.  "Intermediate  elements"  will  be  regarded  as 
geocentric  position  and  velocity,  or  the  standard  elements  a, 
c,  i,  ft  ,  u>,  M  or  some  other  set  of  parameters  which  are 
normally  obtained  in  first  approximation  methods  using  incom¬ 
plete  observations. 

(b)  The  "representation"  or  calculation  of  what  the  observations  would 
be  if  the  basic  intermediate  elements  were  correct.  This  process, 
which  ends  in  the  determination  of  residuals,  may  involve  additional 
observations  as  well  as  those  upon  which  the  intermediate  elements 
are  based. 

(c)  The  differential  correction.  The  differential  correction  utilizes 
residuals  in  differential  formulae  that  relate  them  to  corrections 
of  the  adopted  elements  ^E..  Finally,  we  transform  the  set  of 
intermediate  elements  into  an  adopted  set  of  "terminal  elements" 

E  =  E.  +  AE.  .  Terminal  elements  may  be  initial  position  and 
velocity,  but  more  often  they  will  be  elliptical  parameters  such 
as  a,  e,  i,  ft,  u>,  M  or  some  other  parameters  especially 
suited  in  a  particular  problem  to  ephemeris  computation  and 
analytical  correction  programs. 

It  is  in  the  first  determination  of  the  intermediate  elements  that  one  finds 
the  greatest  number  of  alternative  procedures  or  methods.  This  variety  is  a 
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result  of  mathematical  complexities  and  of  differing  observational  patterns. 

With  observations  scattered  in  time,  kind,  and  quality,  it  is  not  possible  to 
select  a  single  orbit-determination  method  that  will  be  satisfactory  in  all 
circumstances. 

An  orbit  may  be  determined  from  a  six-dimensional  fix,  i.  e. ,  topocentrie 
position  and  velocity.  The  simple  coordinate  transformation  performed  in  this 
method  makes  it  possible  to  determine  the  orbit  without  approximations,  and 
hence,  no  eorreetions  are  needed  until  additional  observational  data  are 
available. 

Two  or  more  complete  three-dimensional  topoeentric  position  fixes  may 
be  available  with  present  radar  equipment.  Geoeentrie  fixes  may  be  obtained 
immediately  from  the  relation 


jl  -  p  -  R 

and  used  as  intermediate  elements.  The  quantities  £  and  Ji_  are  the  vectors 
directed  to  the  satellite  from  the  observer  and  the  dynamical  eenter,  respectively, 
and  R_  is  the  vector  from  the  observer  to  the  eenter  of  the  earth.  If  three  or 
more  fixes  are  available,  it  may  be  preferable  to  determine  the  position  and 
velocity  at  some  eentral  date  by  the  use  of  the  Herrick- Gibbs  formula^  1  and 
use  these  as  intermediate  elements.  If  more  than  three  observations  are  used, 
a  least-squares  reduction  of  the  velocity  at  the  central  date  is  obtained.  Sinec 
the  vcloeity  is  obtained  from  series  expressions  negleeting  higher  terms, 
representation  is  necessary  to  determine  if  differential  eorreetion  is  necessary. 

A  limitation  of  the  Herriek-Gibbs  program  is  its  requirement  for  three 
radar  fixes  separated  by  not  more  than  one  radian  of  arc  as  seen  from  the 
eenter  of  the  earth.  Observations  spaeed  over  any  are  ean  be  processed  by  the 
"Two- Position  Program"  developed  by  Acronutronie;  the  data  may  be  separated 
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by  several  revolutions  of  the  satellite  in  its  orbit.  Observations  spaced  at 
geocentric  angles  ft  and  180  degrees  apart  produce  singularities  in  the  method, 
as  two-position  vectors  with  either  of  these  particular  singularities  make  it 
impossible  to  ascertain  the  orbit  plane. 

Aeronutronic  has  conducted  much  experimentation  with  the  Herrick- Gibbs 
initial  orbit  programs  to  demonstrate  the  usefulness  of  range  rate  measurements 
in  accurate  determination  of  the  orbital  elements,  particularly  in  the  determi¬ 
nation  of  the  period.  The  observations  were  those  made  by  the  Millstone  radar 
over  a  period  of  about  one  week.  The  elements  were  differentially  corrected 
with  all  the  observations  so  as  to  obtain  a  standard  against  which  elements  out 
of  the  initial  orbit  program  could  be  compared.  The  results  are  good,  con¬ 
sidering  the  short  radar  track  of  one  and  two  minutes.  The  use  of  range  rate 
significantly  improved  the  period  determination  in  all  but  one  test  case.  The 
effect  on  the  rest  of  the  elements  is  not  conclusive  and,  in  any  event,  the  error 
in  the  period  is  orders  of  magnitude  more  significant  than  the  errors  in  any 
of  the  other  elements. 

A  great  deal  of  literature  is  available  on  orbit  determination  from  sets  of 
angles  and  angular  rates.  The  basic  methods  are  frequently  associated  with 
the  names  of  LaPlace,  LaGrangc,  Gauss,  and  Gibbs,  though  there  have  been 
many  modifications  by  subsequent  writers.  The  relatively  poor  observations 
of  earth  satellites  limit  the  usefulness  of  some  of  the  classical  methods. 

In  the  LaPlacian  method,  one  determines  position  and  velocity  at  a  central 
date,  to  be  used  as  "intermediate  elements"  for  correction  purposes  prior  to 
the  determination  of  the  "terminal  elements.  "  The  first  approximation  is  based 
upon  Taylor's  series  expansions  of  the  observed  angular  coordinates  or  direction 
cosines;  the  large  angular  rates  further  limit  its  usefulness  for  satellites.  This 
approach  fails  for  close  earth  satellites  because  the  rapid  topocentric  angular 
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motion  gives  rise  to  prohibitive  truncation  errors.  In  order  to  use  the  method 
successfully  for  such  satellites,  one  either  processes  a  large  number  of  angles 
in  a  least-squares  reduction  or  employs  additional  data,  e.g.  ,  range  rate 
observations. 

The  methods  of  LaGrange,  Gauss,  and  Gibbs,  in  their  simplest  and  most 
effective  forms,  all  make  use  of  the  same  first  approximation.  This  approxi¬ 
mation  is  also  based  on  Taylor's  series  expansions,  but  in  the  dynamical 
coordinates  (e.  g. ,  geocentric  coordinates  for  geocentric  orbits)  rather  than  in 
the  observational  ones.  Thus,  for  geocentric  orbits,  not  only  is  the  angular 
motion  greatly  reduced,  but  also  it  is  possible  to  include  some  of  the  higher 
derivatives  in  the  series  with  the  aid  of  the  dynamical  properties  of  motion.  The 
resulting  first  approximation,  which  yields  three  geocentric  position  vectors, 
may  be  more  successful  in  one  of  the  methods  than  in  the  others  because  of  the 
way  in  which  the  position  vectors  are  used  in  the  determination  of  the  terminal 
elements. 

The  first  approximations  of  LaPlace  and  LaGrange- Gauss- Gibbs  share  an 
indeterminacy  that  occurs  when  the  basic  observations  lie  on  a  great  circle  arc. 
The  indeterminacy  may  be  overcome  with  three  observations  if  the  plane  of  the 
great  circle  docs  not  pass  through  the  dynamical  center.  If  it  does,  one  must 
seek  n  fourth  observation  from  which  the  indeterminacy  is  eliminated. 

It  would  be  impossible  to  go  into  the  details  of  these  orbit  computations 
in  this  short  report.  Reference  1  provides  equations  and  procedures  for  use  of 
the  Herrick-Gibbs  method  and  for  those  of  LaPlace,  LaGrange  and  Gauss. 

[  2l 

Moulton  provides  a  thorough  discussion  of  the  methods  of  LaPlace  and  Gauss. 

[  3] 

Herget's  privately  published  book1  provides  detailed  numerical  procedures 

and  manv  numerical  examples  for  computing  orbits  using  the  methods  of  LaPlace 

[  4] 

and  Gauss.  Watson's  book,  though  published  in  1885,  is  still  a  useful  guide. 
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Many  new  methods  for  computing  earth  satellite  orbits  using  canonical  and  other 
element  sets  have  been  proposed.  The  most  significant  recent  works  are  papers 
by  Brouwer,  Kozai,  Garfinkel,  and  Vinti.  ^  6’ 

Once  the  intermediate  orbital  elements  are  determined  from  the  initial 
orbit  routines,  the  approximate  position  can  be  easily  computed  and  compared  to 
the  observed  position  of  the  satellite.  One  thus  gets  a  residual  between  the  observed 
and  computed  position 

-  a  =  A  a 
0  c 

where  a  may  be  range,  altitude,  or  azimuth  (or  right  ascension  and  declination). 

Since  a  is  a  function  of  the  six  orbital  elements  e  ,  e  ,  ...  e  .  We  can 

12  6 

expand  a  into  a  Taylor  series  of  the  parameters  e  ,  ...  e  .  Thus, 

1  6 


A  a  =  an  -  a 
0  c 


order  terms. 


3a  ,  3a 

9  e,  Aei  +  3e„ 


9  o' 

AV---  +  £T; 


Ae  +  higher 

O 


Neglecting  the  higher  order  terms,  we  have  the  equations  for  the  differential 
correction.  Since  the  partial  derivatives  da/de ^  are  functions  of  the  orbital 
elements  and  can  easily  be  computed  (actually,  the  formulas  are  somewhat  long 
and  complicated),  we  have  an  expression  for  the  error  in  the  observations  due 
to  the  errors  in  the  orbital  elements.  What  we  wish  to  find,  however,  are  the 
errors  in  e, ,  e0,  ...  e  due  to  the  differences  in  the  observed  and  computed 
positions.  Therefore,  if  we  have  six  observed  differences  Aa  ,  Ao'g’  ••• 
we  can  invert  the  matrix  of  coefficients  of  the  Ae.'s: 
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9o^  9g^  9a^ 

Ac*  =  - —  Ae,  +  - —  Ae  +  . .  .  + -  Ac 

19c  1  9c  2  de„  6 

12  6 


da2  9«2 

A  a.  =  - — -  Ae,  +  - —  A  e  +  . . .  +  -  A  e„ 

2  8  1  9e2  2  9eo  G 


8"g  8"g  3“g 

a“g  =  5T  Aei  V  AV  •••  Aeo 

12  G 


and  solve  for  the  corrections  to  the  elements,  the  Ae.  . 

If  wc  have  more  than  six  observational  residuals  and  hence,  more  than 
six  equations  of  condition,  then  we  can  perform  a  least-squares  fit.  In  matrix 
form,  we  have 

A  X  =  R  , 

where 
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and  m  >  n  . 

Usually  n  =  6  (the  number  of  orbital  elements),  but  sometimes  a  seventh 
"element"  or  drag  parameter  is  used  in  the  differential  correction. 

It  is  appropriate  to  point  out  that  each  equation  of  condition 


da  A  3  a  A  3a. 

Aa  =  - —  Ae,  +  - —  Ae  +  . . .  +  - —  Ae„ 
3e  1  3e  2  3e„  6 

1  Z  o 
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may  be  weighted  by  multiplying  it  by  1/ct  where  c r  is  the  root-mean-square 
(r.  m.s. )  of  the  error  of  the  observation.  The  values  of  cr  are  different  for 
different  radars  and  various  types  of  optical  observations  and  are  determined 
independently  of  the  differential  correction  process.  Quite  often,  however, 
since  the  weighting  of  these  equations  does  not  contribute  substantially  to  the 
accuracy  of  the  final  results,  and  since  the  values  of  a  are  often  difficult  to 
determine,  a  common  practice  is  to  omit  the  weighting  altogether  and  assign 
each  equation  a  unit  weight  of  one  (e.  g. ,  SPADATS). 

Since  little  is  known  about  the  values  of  a  for  the  various  tracking 
installations,  a  study  of  this  problem  is  being  conducted  at  Aeronutronie.  NASA, 
on  the  other  hand,  does  make  use  of  a  weighted  least-squares  technique  in  the 
Project  Mercury  Program. 

In  matrix  form,  we  may  get  the  normal  equations  (usually  six)  by  multi¬ 
plying  by  the  transpose  of  A. 

T  T 

A  A  X  =  A  R. 

B  X  =  C 

We  will  then  have  six  equations  in  six  unknowns  which  can  be  solved  by  standard 
elimination  methods  such  as  the  Gauss- Jordan  reduction.  The  solution  of  these 
equations  will  give  us  the  correction  to  the  elements  Ae^,  Ae^,  .  . .  Ae^.  . 

At  Space  Track,  after  forming  a  set  of  residuals  corresponding  to  a  set 
of  observations,  the  magnitudes  of  the  angle  and  range  residuals  are  compared 
with  an  absolute  maximum  value  of  1000  km,  and  the  range  rate  residual  is 
compared  with  an  absolute  maximum  value  of  0.5  km/sec.  All  residuals 
exceeding  these  limits  are  rejected  from  the  current  iteration.  Upon  completion 
of  the  first  residual  rejection  test,  a  second  similar  residual  rejection  test  is 
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performed  with  a  new  maximum  value  equal  to  (1. 5)  (r.  m.  s.  value  of  the 
previously  accepted  residuals).  The  remaining  number  of  residuals  is  then 
used  in  the  differential  correction  equations.  Upon  completion  of  the  calculation 
of  the  new  set  of  corrected  elements,  however,  the  absolute  maximum  value  of 
the  first  rejection  test  is  subject  to  change  from  1000  to  75  km  if  the  new  r.  m.  s. 
value  of  the  residuals  is  less  than  50  km.  With  the  existing  system  of  not  too 
accurate  radars,  it  is  necessary  to  have  many  observations  to  get  a  good  satellite 
orbit.  It  is  necessary  to  use  a  large  number  of  observations  (usually  20  or  more 
per  satellite)  in  a  "least-squares  differential  correction"  every  one  or  two  weeks 
to  update  a  satellite  orbit;  however,  satellite  orbits  not  appreciably  affected  by 
drag  can  often  be  accurately  predicted  without  running  the  differential  corrections 
for  periods  of  from  30  to  60  days. 

The  literature  contains  three  excellent  and  most  complete  references  on 

[  4  9  10] 

the  differential  correction  methods.  ’  ’  References  4  and  10  contain 

excellent  sections  on  the  least-squares  technique  and  on  the  weighting  of  residuals, 
though  the  notation  is  somewhat  clumsy  since  matrix  notation  is  not  used. 

One  should  point  out  that  the  particular  set  of  orbital  parameters 

ei’  e2’  ' ' '  ’  e6  °nC  c^ooses  anc*  the  particular  method  of  orbit  computation 

depend  on  the  kind  of  orbit.  Since  the  perigee  cu  is  not  well-determined  in  a 

nearly  circular  orbit,  it  is  necessary  to  use  other  orbital  element  seta  which 

can  be  better  determined.  For  orbits  of  low  eccentricity,  the  elements  a, 

e  cosu,  esinw,  ft,  i,  M  are  useful.  Since  the  node  ft  is  not  well- 

determined  in  a  nearly  equatorial  orbit,  it  is  also  necessary  to  use  a  different 

orbital  element  set  in  this  case.  For  orbits  of  low  inclination,  the  elements 

a,  e,  u>,  sin  i  sin  ft,  sin  i  cos  ft,  M  are  sometimes  used.  For  orbits  of  low 

inclination  and  low  eccentricity,  Herrick  has  found  orbital  elements  that  depend 

on  vector  quantities  such  as  U  ,  V  most  useful,  where  U  is  a  unit  vector 

— o  — o  — o 
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in  the  direction  of  the  satellite  at  time  t  and  V  is  perpendicular  to  U 

o  —  o  —  o 

and  in  the  orbital  plane  of  motion  of  the  satellite.  A  differential  correction 
developed  Eckert  and  Brouwer^  11  ^  makes  use  of  infinitesimal  rotations  and 
corrects  on  A  tp  ,  A  ip  ,  and  Aip  where  these  are  the  resultants  of  rotations 

1  Z  d 

of  the  orbit  about  the  x,  y,  and  z  axes,  respectively.  Actually  the  Eckert- 
Brouwer  method  uses  combinations  of  these  and  other  elements  to  yield  several 
sets  of  elements,  each  suited  to  a  particular  type  of  orbit  (for  example  AM  + 
Aip  ,  Aip  ,  Aip  ,  eA ip  ,  Aa/a,  and  Ae  are  used  for  low  eccentricity  orbits). 

u  1  Z  U 

All  in  all,  there  arc  many  different  sets  of  orbital  elements  and  computation 
methods,  each  with  its  own  advantages  and  disadvantages. 
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SECTION  IV 

PERTURBA1  IONS  OF  THE  SATELLITE  ORBIT 

The  significant  perturbative  forces  are  oblateness  of  the  earth,  drag,  and 
solar  radiation  pressure.  Less  significant  perturbative  forces  are  those  due 
to  magnetic  effects,  lunar  and  solai  gravitational  attraction. 

Oblateness 

The  potential  function  for  the  earth  can  be  written  as 


2  .  2 . .  3  .  3  _  .  „ . 

— —  (3  sin  6-1) - —  (5  sin  o  -  3  sin  6) 

2r  2r 

J4  4  2 

- —  (35  sin  6  -  30  sin  6  +  3) 

8r 

*^5  5  3 

— Z  i63  sin  6-70  sin  6+15  sin  6) 

8r° 

6-  (231  sin66  -  315  sin46  +  105  sin~6  -  5  )  > 

16r6  J 


where 

6  is  the  instantaneous  latitude  of  the  satellite, 
p  is  the  gravitational  constant  for  the  earth  1.407639, 
r  is  the  radial  distance  of  the  satellite  from  the  center  of  the 
earth  in  earth  equatorial  radii, 

J9  =  1.08228xl0-3  (±  .  00003xl0-3), 

-6  ct 

J  =  -2.3x10  (±.2xl0-b), 

=  -2. 12xl0-6  (±  .  05x1 0-6), 
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J„  =  -.2x10_G(± . lxlO'6),  and 

5 

J„  =  1.0xl0’G  (±.8x10'G). 

6 


For  most  practical  estimates  of  the  oblateness  effeet,  the  J>}  term  is 
the  only  one  that  need  be  considered.  These  oblateness  effeets  are  quite  con¬ 
siderable  on  close  earth  satellites  and  ean  affeet  the  perigee  and  nodal  points 
by  as  much  as  several  degrees  per  day.  One  ean  see  by  looking  at  the  potential 
function  that  as  r  increases,  the  effect  of  the  oblateness  is  diminished. 

The  oblateness  of  the  earth  primarily  affects  three  orbital  elements  £2, 

w,  and  M;  however,  the  small  effeet  on  the  remaining  elements  ean  be 

[  5  6  7  8  1 2 1 

accurately  determined  analytically.  ’  ’  ’  ’  The  equatorial  bulge  eauses 
the  node  £2  to  preeess  in  a  manner  similar  to  that  of  a  spinning  top. 

The  formula  for  this  preeession  is  given  by 

J'  n  eos  i 

n=-  — - 

a2(l  -  e2)2 

where 

J’  =  1. 62345x10" 3 , 

a  =  the  semi-major  axis, 
e  =  the  eccentricity, 

i  =  the  inclination  of  the  orbit  to  the  equator,  and 
n  =  the  mean  motion  (angular  velocity)  of  the  satellite 
in  its  orbit. 


We  should  note,  however,  that  at  i  =  90  degrees,  cos  i  =  0  and  hence, 

£2=0;  also,  £2  can  reach  a  maximum  of  slightly  more  than  10  deg/day  for 
certain  low  inclination  orbits.  Actually,  £2  is  referred  to  as  the  regression 
of  the  node  since  its  sign  is  minus.  Figures  3  and  4  are  most  instructive. 
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INCLINATION  ANGLE,  I  ( deg) 


Fig .  4  Nadal  Regression  Rate  as  Functions  of  Mean  Altitude  and 
Orbital  Inclination  Angle  (from  Reference  13) 
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The  second  important  effeet  of  the  bulge  at  the  earth's  equator  is  that  the 
orbital  ellipse  rotates  steadily  in  its  own  plane  (see  Fig.  5).  The  rotation  rate 
is  given  by  the  formula 

J  '  n  (5  eos“  i  -  1 ) 

•  2 _ 


The  rotation  is  forward,  in  the  same  direction  as  the  satellite  motion,  for 
i  <  63.4  degrees,  zero  when  i  =  63.4  degrees,  and  backward  for  i  >  63.4 
degrees.  We  should  note  from  Fig.  6  that  this  preeession  of  perigee  (often 
called  apsidal  precession)  can  reach  20  degrees  in  certain  near  equatorial 
orbits. 

ORBITAL 


Fig.  5  Diagram  Showing  Rotation  of  the  Orbital 
Ellipse  in  Its  Own  Plane.  (Initially  the 
orbit  is  nearest  to  the  earth  at  A,  and 
later,  at  B.) 
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Fig  .  6  Apsidal  Precession  Rate  as  a  Function  of  Mean  Altitude 
and  Orbital  Inclination  Angle  (from  Reference  13) 
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The  mean  anomaly  M  of  the  orbit  is  also  severely  perturbed  by  the 
earth's  oblateness.  Thus, 


M  = 


J  ^n  (3  cos  i  -  1 ) 
2a2  (l-e2)f 


The  oblateness  has  a  small  effect  on  the  remaining  orbital  elements,  but  they 
are  generally  not  too  important  for  short  time  periods.  Complete  investigations 
of  the  effect  of  oblateness  on  all  the  orbital  elements  may  be  obtained  from  the 
papers  of  Brouwer,  Kozai,  and  Garfinkel.^’  ^ 


Drag 

The  perturbative  effects  of  air  drag  are  significant  at  altitudes  below 

about  300  miles.  It  is  usual  to  express  the  aerodynamic  drag  D  of  a  satellite 

in  terms  of  a  drag  coefficient  ,  based  on  the  maximum  frontal  area  A 

D  c 

of  the  satellite: 


D  =  —  C  A  pV2 
2  D  cH 


where 

p  is  the  local  air  density  and 
V  is  the  magnitude  of  the  satellite  velocity. 


In  most  formulations,  however,  the  drag  deceleration  is  found  to  be  the  basic 
quantity 


D  _  CDAe  2 
m  2m  P 


2 

Bp  v  • 
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Considerable  simplification  may  be  achieved  by  the  use  of  the  ballistic  coefficient 
B,  since  the  problem  may  be  solved  parametrically  for  a  series  of  values  of 
B  which  eliminates  the  necessity  for  variables  C  ,  A^,  and  m. 

For  a  given  vehicle  mass  m  is  obviously  known,  but  for  more  complicated 
shapes  each  particular  case  will  necessarily  determine  its  own  coefficient  of 
drag  and  reference  area.  The  drag  coefficient  for  upper  atmospheres  is  a  func¬ 
tion  of  the  flow  regime  and  the  body  shape  wherever  the  reference  area  is  a 
function  of  shape  and  body  dynamics  such  as  orientation  tumbling,  etc.  It  can  be 
shown  that  at  about  75  statute  miles  the  satellite  is  out  of  the  slip  flow  region, 
and  free  molecular  flow  conditions  apply  for  any  velocity  and  any  satellite  size. 
For  free  molecular  flow,  the  fraction  of  diffusely  reflected  particles  is  0.9; 
thus,  the  curve  for  diffuse  reflection  should  be  used  for  design  consideration. 

It  is  noted  in  Reference  13  that  for  high  molecular  speed  ratios,  a  2  is 
reached  for  both  spherical  and  cylindrical  satellites.  This  value  is  a  valid 

first  approximation.  The  reference  area  for  steady  flight  is  the  cross-sectional 

7T  2 

area  perpendicular  to  the  air  stream,  i.e.  ,  for  a  sphere  =  —  D  .  For 
tumbling  satellites  it  is  assumed  that  each  orientation  is  equally  probable,  and 
the  average  reference  area  is  given  by  the  total  surface  area  divided  by  four. 
Numerical  results  from  Reference  13  arc  compared  in  Fig.  7  for  two  satellites. 

It  can  be  seen  that  for  a  cylinder  with  a  diameter-to-length  ratio  of  0.  078,  the 
relative  drag  for  random  tumbling  is  11.9  times  higher  than  for  an  orientation 
with  the  main  axis  parallel  to  the  air  stream,  and  28  percent  lower  than  for  an 
orientation  with  the  axis  perpendicular  to  the  air  stream. 

In  addition  to  the  satellite  drag  based  strictly  on  the  air  densities  of  an 
atmospheric  model,  such  as  ARDC  1959,  certain  special  effects  exist  which 
should  be  considered. 
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Aeaumed  Surface  Interaction:  Dlffuae  Re-ernlealon 

Accommodation  Coefficient  -  0.50 
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A  Approximated  by  drag  coefficient  of  cylinder  with  axla  normal  to  stream. 

'Normalized  on  the  baee  area  Jd^  -  0.212  ft. 

Fig.  7  Effect  of  Attitude  and  Tumbling  on  Drag  of 

Representative  Sottellites  (from  Reference  13). 


The  earth's  oblateness  basically  adds  a  latitude  effect  into  the  atmospheric 
tables,  which  is  presently  missing  in  all  the  existing  air  density  models.  Good 
discussions  of  the  oblateness  effects  on  the  air  drag  are  given  in  the  Literature, 
References  14,  15  and  16. 

Rotation  of  the  atmosphere  causes  the  orbital  lifetimes  of  satellites  launched 
in  retrograde  orbits  to  be  different  from  those  of  satellites  in  corresponding  direct 

f  j 7  1 

orbits.  This  problem  has  received  some  attention.  ’ 


23 


TM-3641 


Diurnal  atmospheric  bulge  exists  above  125  statute  miles  and  is  believed 
to  result  from  thermal  and  gravitational  air  expansion  effeets.  Also,  the 
atmosphere  is  very  unstable  between  50  and  70  statute  miles  due  to  large  waves 
and  tidal  motions  in  the  air.  At  125  statute  miles,  the  solar  flares  may  cause 
ehanges  of  40  pereent  in  air  density,  becoming  more  important  above  175 
statute  miles.  There  seems  to  be  large  variations  related  to  the  seven-year 
sun  spot  eyele,  whieh  would  imply  that  the  satellite  lifetime  is  also  a  function 
of  the  year  in  whieh  the  orbit  is  established. 

Sinee  it  is  well-known  that  drag  is  the  major  influence  on  satellite  life¬ 
time,  it  is  appropriate  to  mention  that  a  fairly  accurate  estimate  of  a  satellite's 
lifetime  ean  be  obtained  at  the  beginning  of  its  eareer  if  the  initial  eccentricity 
e^  is  known  and  is  not  larger  than  about  0.  2,  and  if  the  period  of  revolution  has 
been  measured  for  several  days.  The  total  lifetime  of  the  satellite,  t  ,  is 

■Lj 

given  by 


where 

T  is  the  initial  period  of  revolution,  and 
o 

x  is  the  daily  decrease  in  the  period  of  revolution. 

If,  for  example,  e^  =  0.1,  Tq  =  100  minutes  and  x  =  0.  05  minutes/day, 

we  find  t  =  150  days,  and  this  estimate  of  lifetime  should  not  be  in  error  by 
L 

more  than  about  15  days. 

The  most  important  faets  to  remember  about  air  drag  are  that  it  varies 
widely  and  it  is  difficult  to  prediet  these  short  period  fluctuations.  Figures 
8,  9  and  10,  dealing  with  air  density  and  drag,  should  prove  useful  and 
informative. 
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Fig.  8  Chart  of  the  Atmosphere  Below  200  Miles  Altitude,  Showing  Air 
Density,  Temperature,  and  Concentrations  of  Molecules  and  of 
Free  Electrons  (from  Reference  19) 
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The  points  plotted  are  the  values  of  air  density  obtained  from  the 
orbits  of  each  of  the  satellites,  with  the  full  line  drawn  through 
them  representing  the  average  density.  The  broken  line  shows  the 
density  according  to  the  pre-satellite  ARDC  model,  which  is  too 
low  by  a  factor  of  between  3  and  1  1  at  heights  between  120  and 
250  miles.  Sea-level  density  is  0.0765  Ib/cu.  ft.,  or  0.00123 

grams/cm^. 


Fig.  9  The  Variation  of  Air  Density  with  Height  (from  Reference  19) 
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The  Figure  shows  that  density  departs  from  its  average  value  by  up  to 
35  percent,  and  that  the  maximum  values  of  density  occur  at  intervals 
of  about  28  days. 


Fig.  10  The  Variation  of  Air  Density  during  1958,  at  Heights 
between  100  and  150  Miles,  as  Given  by  Sputnik  3 
Rocket  (from  Reference  19) 
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Aerodynamic  drag  primarily  affects  the  elements  M,  a,  and  e.  In  a 
circular  orbit,  the  air  drag  slows  down  the  speed  of  a  satellite  slightly,  but  as 
its  speed  drops  below  the  orbital  speed  proper  to  its  height,  it  begins  to  descend 
at  a  small  angle.  As  soon  as  this  happens,  gravity  aeeelerates  it  on  its  down¬ 
hill  path.  So  the  satellite  begins  to  go  faster,  and  the  speed  increases  beyond 
its  original  value  until  it  reaches  the  slightly  higher  speed,  which  is  the  proper 
orbital  s ">eed  at  this  new  and  slightly  lower  height.  Then  drag  begins  to  slow 
it,  and  the  whole  cycle  begins  again.  The  net  result  is  that  the  satellite  slowly 
descends  at  a  steadily  increasing  speed  in  a  nearly  circular  spiral.  In  this  case, 
the  semi-major  axis  a  is  continually  decreasing  and  since  the  period 


the  period  also  decreases. 

If  the  initial  orbit,  instead  of  being  circular,  is  appreciably  elliptic,  the 
effect  of  aerodynamic  drag  is  quite  different  at  first.  Since  air  density  falls  off 
rapidly  as  height  increases,  a  satellite  in  an  elliptic  orbit  will  get  a  much  greater 
wallop  from  drag  at  or  near  perigee  than  at  or  near  apogee.  In  faet,  the  drag 
can  be  largely  ignored  exeept  in  the  region  near  perigee.  Over  this  short  seetion 
of  the  orbit,  drag  causes  a  small  loss  in  speed  of  the  satellite.  The  minimum 
height  of  the  satellite  is  reduced  only  to  a  very  negligible  degree,  but  the  maxi¬ 
mum  height  at  apogee  is  reduced  to  a  much  greater  extent.  The  effect  of 
aerodynamie  drag  on  an  elliptic  orbit,  therefore,  is  to  make  the  orbit  more 
nearly  circular  by  steadily  reducing  the  maximum  height  and  scarcely  reducing 
the  minimum  height.  As  in  the  circular  case,  the  semi-major  axis  a  and  the 
period  P  are  also  decreasing.  In  eccentric  orbits,  however,  the  eccentricity 
is  also  decreased  since  the  orbit  continually  grows  more  nearly  eireular. 
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Very  gross  errors  in  satellite  position  and  the  orbital  elements  would 
develop  in  a  few  orbits  if  the  major  perturbations  sueh  as  oblateness  and  drag 
were  ignored.  One  ean  fairly  easily  account  for  the  major  changes  in  earth 
satellite  orbits  eaused  by  oblatencss  by  using  the  formulas  for  the  seeular  vari¬ 
ations  in  fi.w,  and  M  as  given  in  the  seetion  on  oblateness.  The  effeet  on 
the  mean  anomaly  M  by  drag  is  a  very  considerable  one.  If  one  writes  the 
mean  anomaly  as  a  function  of  time 

2  3 

M  =  Mq  +  c^  (t  -  T)  +  e  (t -  T)  +  Cg  (t  -  T)  +  periodic  terms  . 

then  the  e^  constant  is  just  the  Kepler  two-body  change  in  the  orbit  plus  the 
secular  perturbations  due  to  oblateness.  The  e  and  c  constants  which  are 

L  O 

the  drag  terms  are  usually  determined  empirically  by  a  least-squares  fit,  i.  e. , 
the  differential  corrections  corrects  on  the  seven  elements  a,  e,  i,  SI  ,  cu, 

M,  and  c  and  e  are  only  used  in  the  ease  of  very  high  drag  satellites; 

^  O 

Smithsonian  and  SPASUR  have  found  that  the  use  of  the  c  term  with  low  drag 

o 

satellites  often  decreases  rather  than  increases  the  accuracy  of  prediction. 

Drag  correeiions  based  on  this  empirieal  c^  can  also  be  made  to  the  elements 
a  and  e  in  a  straightforward  manner: 

a  =  ao  '  1  aoC2(t"T)  ’ 

and 

a  (1  -  e  ) 

,  o  o 

e  =  1  -  - 

a 

All  existing  operational  systems  in  NASA,  Space  Track,  NAVSPASUR, 
JPL,  Smithsonian,  ete. ,  make  use  of  these  simple  corrections  for  drag  and 
oblatencss;  however,  the  inelusion  of  more  terms  (i.e. ,  short  and  long  periodic 
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[  5] 

perturbations  and  more  sophisticated  techniques  is  quite  prevalent  in  the 
above  tracking  and  computational  centers.  Figures  11  through  17  serve  to 
illustrate  and  clarify  this  discussion. 


MINIMUM  HEIGHT  IS  REDUCED 
ONLY  VERY  SLOWLY 


Fig .  1 1  Shrinking  of  Satellite  Orbit  Under  the  Action  of  Air 
Drag  —  Not  Exactly  to  Scale  (from  Reference  19) 


Solar  Radiation  Pressure 

Solar  radiation  pressure  produces  significant  perturbations  only  on  large 

light  satellites  such  as  Echo.  More  specifically,  this  perturbation  is  significant 

2  [  13 

for  satellites  with  area-to-mass  ratios  of  approximately  25  em  /gm  or  greater,1 
although  solar  radiation  pressure  also  produces  measurable  effects  on  small, 
dense  satellites  at  high  altitudes.  Solar  radiation  pressure  is  estimated  to  be  equal 
to  the  air  drag  at  a  500-statute  mile  altitude.  Its  must  important  effect  is  dis¬ 
placing  the  geoeentric  eenter  of  the  orbit.  For  eertain  resonanee  conditions,  the 
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TIME  FROM  RERIBEE  -  MINUTES 


The  average  drag  during  the  20  minutes  when  drag  is 
appreciable  is  just  under  50  percent  of  the  maximum 
drag,  the  area  af  the  rectangle  shown  being  equal  to 
the  area  under  the  curve.  The  graph  applies  to  a 
satellite  with  orbital  eccentricity  near  0. 1 . 


Fig.  12  The  Variation  in  the  Drag  Acting  on  a  Satellite  Near 
Perigee  (fram  Reference  19) 
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Numbers  on  the  curves  indicate  the  initial  eccentricity  e  . 


Fig.  13  Variation  af  the  Eccentricity  e  of  the  Orbit  during  a 
Satellite's  Lifetime  (from  Reference  19) 
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Fig.  14  Altitude  Decay  History  for  a  Circular  Orbit  (from  Reference  13) 
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Fig  15  Lifetimes  for  Circular  Orbits  (from  Reference  13) 
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INITIAL  PERIGEE  ALT  ITU  DE  (STATUTE  MILES) 
Fig.  16  Satellite  Lifetimes  (from  Reference  13) 
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17  Satellite  Lifetimes  as  Functions  of  Apogee  and  Perigee 
Altitudes  (from  Reference  13) 


TM-3641 


orbital  lifetime  is  changed  by  significant  amounts.  For  the  Beacon  satellite 

2 

(area-to-mass  ratio  of  23.  2  cm  /gm)  at  a  mean  altitude  of  950  miles,  the  life¬ 
time  can  vary  by  a  factor  of  10,  depending  on  the  hour  of  launch. 

In  general,  during  a  complete  orbital  period,  solar  radiation  pressure 
causes  a  first-order  perturbation  of  all  six  orbital  parameters.  However,  the 
most  conspicuous  effect  for  a  nearly  circular  orbit  is  a  displacement  of  its 
geometric  center.  This  displacement  is  perpendicular  to  the  earth-sun  line  in 
the  orbit  plane  and  in  a  direction  such  as  to  decrease  the  altitude  of  that  part 
of  the  orbit  in  which  the  satellite  moves  away  from  the  sun.  Calculations  show 
that  at  a  mean  altitude  of  1000  miles,  radiation  pressure  can  displace  the  orbit 
of  the  100-foot  Echo  balloon  at  rates  up  to  3.  7  miles/day,  the  orbit  of  the  12- 
foot  Beacon  satellite  at  0.  7  milcs/day,  and  even  the  orbit  of  Vanguard  I  at  a  much 
slower  rate  of  about  one  mile/year. 

For  certain  resonance  conditions,  these  perturbations  due  to  solar  radi¬ 
ation  pressure  accumulate  and  drastically  affect  the  satellite's  lifetime.  For 
the  Beacon  satellite  at  a  mean  altitude  of  950  miles,  an  initial  eccentricity  of 
0. 106  and  an  inclination  of  40  degrees,  the  lifetime  can  vary  by  a  factor  of  10, 
depending  on  the  hour  of  launch.  For  an  inclination  of  48  degrees  of  the  Beacon 
satellite,  these  conditions  are  no  longer  resonant,  and  the  variation  in  lifetime 
is  reduced  to  a  factor  of  2. 

For  the  Echo  balloon  placed  in  a  1000-milc  altitude  circular  orbit  and  an 
inclination  of  35  degrees,  we  find  that  the  lifetime  is  240  days.  For  an  initially 
circular  equatorial  orbit,  the  resonance  altitude  of  4000  miles  leads  to  a  1.3- 
year  lifetime,  while  the  same  orbit  at  1000  miles  altitude  has  an  extremely 
long  lifetime. 

For  an  excellent  report  on  the  effects  of  solar  radiation  pressure  on 

[20  21 1 

satellite  orbits,  one  is  urged  to  read  papers  of  1. 1.  Shapiro.  ’  They  are 
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brief,  lucid  and  give  very  valuable  information.  For  the  present,  we  offer 
Figures  18  through  23,  which  depict  the  radiation  effects  on  eccentricity, 
perigee  height,  and  the  argument  of  perigee,  a.' . 


DAYS  SINCE  VERNAL  EQUINOX - ► 


Fig.  18  Time  Variation  of  Perigee  Height  ond  Meon  Altitude 
for  Echo  I  (from  Reference  21) 
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Fig.  19  lime  Variation  of  Argument  of  Perigee 
(from  Reference  21) 


SHADOW 


TM- 3(541 


Fig.  22  Time  Variation  af  Perigee  Height  and  Mean  Altitude 
far  the  Beacon  Satellite  (from  Reference  20) 


The  data  are  published  orbital  elements  with  effect 
af  the  third  harmonic  removed. 


Fig.  23  Perturbations  an  Perigee  Height  of  Vanguard  I 
(from  Reference  23) 
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Magnetic  Effeets 

If  the  outer  shell  of  the  satellite  becomes  eleetrieally  charged,  its  move¬ 
ment  in  a  plasma  creates  three  kinds  of  electromagnetic  drag  effeets:  eouloumb 

[  131 

drag,  induction  drag,  and  wave  drag.  It  has  been  estimated1  J  that  for  a 

6  3 

negative  potential  of  10  volts  and  a  partiele  density  of  10  eleetrons/cm  ,  the 
ratio  of  eleetrieal-to-aerodynamie  drag  may  be  about  0.4  for  eertain  eases, 
while  induetive  and  wave  drag  appear  to  be  negligible.  At  altitudes  above  about 
750  statute  miles,  the  magnetic  drag  may  exeeed  neutral  drag  for  large  balloon 
satellites.  Generally,  this  souree  of  drag,  however,  has  a  negligible  effeet. 

Lunar  and  Solar  Perturbations 


Lunar  and  solar  attractions  are  the  major  sourees  of  perturbations  above 
the  24-hour  orbit  altitude  of  about  22,300  statute  miles,  but  their  effects  are 
almost  negligible  for  close  earth  satellites.  The  principal  effect  is  a  regression 
of  the  satellite  orbit  plane  about  the  normal  to  the  orbit  plane  of  the  perturbing 
body.  This  regression  is  given  by 
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where  the  subscript  d  indicates  a  parameter  of  the  disturbing  body.  This 
effect  is  plotted  in  Figs.  24  and  25  for  the  sun  and  the  moon.  There  is  also  a 
radial  or  tidal  perturbation,  the  maximum  value  being  about 


Ar 


max 


}td 


for  circular  orbits  of  radius  r  .  This  value  is  plotted  in  Fig.  26.  A  summary 
graph  of  the  major  perturbations  is  given  in  Fig.  27. 


A/*  C'-* 

R.  H.  Greene 
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Fig.  24  Average  Regression  Rate  about  the  Normal  to  the  Lunar  Orbit 
Due  to  Lunar  Attraction  (from  Reference  13) 
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Fig .  25 


INCLINATION  TO  EQUATOR  td«g) 

Average  Ecliptic  Nodai  Regression  Rate  Due  ta  Solar 
\ttraction  (from  Reference  1 3) 
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Fig  26  Maximum  Radial  Perturbation  Due  to  Attraction  of  the 
Sun  and  Maan  (from  Reference  13) 
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